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Resumen 


Para estimar la disponibilidad del agua subterránea y su evolución es indispensable 
medir los niveles del agua subterránea por periodos largos de tiempo. Sin embargo, 
los recursos económicos son una limitante en la recopilación de suficientes datos 
que permitan estimar con certeza la evolución de dichos niveles. Por este motivo 
se requiere de métodos para el diseño de redes de monitoreo de los niveles del 
agua subterránea que faciliten utilizar los recursos económicos disponibles para 
el monitoreo de manera efectiva. El método de Herrera y Pinder (Herrera, 1998) 
para el diseño óptimo de redes de monitoreo de la calidad del agua subterránea, 
aplicado por Herrera et al. (2001), y también por Herrera y Pinder (2005), tiene 
varias ventajas. Por un lado, permite seleccionar tanto posiciones como tiempos 
de monitoreo en forma óptima; además, con la información obtenida de la red de 
monitoreo, permite hacer estimaciones en espacio y tiempo, y utiliza un modelo de 
flujo y transporte, el cual sintetiza el conocimiento de las características del acuífero, 
de la dinámica del flujo subterráneo y de la problemática de calidad del agua 
subterránea que se tenga, así como del conocimiento de las leyes de la física. En 
este trabajo, una modificación de esta metodología se prueba por primera vez para 
el diseño de redes de monitoreo de los niveles del agua subterránea. Se presenta un 
caso sintético en el que la red se diseña para un solo tiempo, suponiendo que el flujo 
está en estado estacionario. Los resultados presentados indican que la metodología 
propuesta por Herrera y Pinder para el diseño de redes de monitoreo tiene buenas 
perspectivas para aplicarse exitosamente al diseño de redes de monitoreo de los 
niveles del agua. 
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Introducción subterránea y de las cargas de contaminantes. 


Sin embargo, se requiere de un alto costo para 


El diseño de redes de monitoreo del agua 
subterránea consiste en establecer un programa 
de supervisión continua en pozos de muestreo 
para medir niveles del agua subterránea y tomar 
muestras para el análisis de variables físicas, 
químicas y biológicas previamente establecidas. 
Este tipo de redes de monitoreo piezométricas 
y de calidad del agua subterránea son indis- 
pensables para conocer el comportamiento de 
un acuífero, así como sus tendencias de calidad 
del agua. Esta información es básica para lograr 
una gestión eficaz y controlar los impactos 
ambientales negativos de la extracción del agua 


medir los niveles del agua subterránea, recabar 
las muestras, y llevar a cabo el análisis de los 
parámetros de la calidad del agua subterránea 
en todos los pozos de monitoreo de un acuífero 
o de alguna zona contaminada. Es por ello que 
se requiere de un diseño óptimo de redes de 
monitoreo con las cuales se optimice la selección 
de los pozos y los tiempos de muestreo (en 
dónde, cuándo y cuántos pozos muestrear), los 
cuales proporcionen la mayor información para 
lograr estimaciones confiables de los parámetros 
de calidad del agua subterránea y de los niveles 
del agua subterránea. 
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Como respuesta a estas necesidades, He- 
rrera (1998) propuso un método para el 
diseño óptimo de redes de monitoreo, que 
involucra en el diseño el espacio y el tiempo 
de forma combinada. Este método tiene varias 
ventajas: por un lado permite seleccionar 
tanto posiciones como tiempos de monitoreo 
en forma óptima; con la información ob- 
tenida de la red de monitoreo permite 
hacer estimaciones en espacio y el tiempo; 
y utiliza un modelo de flujo y transporte, 
el cual sintetiza el conocimiento de las ca- 
racterísticas del acuífero, de la dinámica del 
flujo subterráneo y de la problemática de 
calidad del agua subterránea que se tenga, a 
la vez que incluye el conocimiento de las leyes 
físicas. Este método se probó en un caso de 
estudio que se reporta en Herrera et al. (2001), 
y otros ejemplos de aplicación se presentan 
en Herrera y Pinder (2005). Nos referiremos 
a este método como el método de Herrera y 
Pinder. 

El método de Herrera y Pinder se pro- 
puso originalmente para diseñar redes de 
monitoreo de calidad del agua y se divide 
en dos procedimientos: 1) estima la varianza 
del error de la estimación que se obtendría si 
se tomaran muestras de agua en diferentes 
pozos y tiempos de monitoreo, y 2) escoge las 
posiciones y tiempos de monitoreo que mi- 
nimizan una función de la varianza del error 
de la estimación, con los cuales se definen la 
red de monitoreo y su programa de muestreo. 
En este trabajo se prueba una modificación de 
ese método para diseñar redes de monitoreo 
de los niveles del agua subterránea. Dicha 
modificación consiste en emplear la carga 
hidráulica en vez de la concentración, y por 
tanto un modelo de flujo en vez de un modelo 
de transporte, como en los dos procedimien- 
tos descritos anteriormente. La formulación 
del método resultante se presenta para redes 
de monitoreo espacio-temporales y la apli- 
cación se realiza en un caso de estudio sintético 
en el que la red de monitoreo se diseña para 
un solo tiempo. 


Estado del arte en las redes de monitoreo 
piezométrico 


A continuación se presenta un resumen de 
trabajos que proponen métodos para el diseño 
de redes de monitoreo del agua subterránea. 
Se incluyen también algunos que presentan 
métodos para el diseño de redes de monitoreo 
de la calidad del agua subterránea, ya que 
con algunas modificaciones pueden servir de 
base en el diseño de redes que permitan una 
adecuada caracterización de la carga hidráu- 
lica en acuíferos. 

Tres puntos de vista han influenciado 
mucho el diseño de redes de monitoreo 
del agua subterránea (ASCE, 2003). En el 
primero de ellos, que llamamos marco hidro- 
lógico, siguiendo a Loaiciga et al. (1992), la 
red y su programa de muestreo se definen 
por consideraciones que toman en cuenta 
solamente las condiciones hidrológicas del 
sitio, sin recurrir a técnicas estadísticas O 
probabilísticas avanzadas. El segundo propone 
el análisis de datos dentro de un marco 
estadístico y define la red de muestreo basado 
en inferencias obtenidas de los datos, éste se 
llama marco estadístico. Al último punto de 
vista le llamaremos aquí marco de modelación, 
en donde se utilizan modelos matemáticos del 
agua subterránea para determinar posiciones y 
frecuencias de muestreo. 

A continuación, primero se presentan los 
trabajos basados en criterios geoestadísticos, 
que pueden ser incluidos en el marco estadís- 
tico, y posteriormente trabajos que pueden 
incluirse en el marco de modelación. En cada 
una de estas categorías primero se describen los 
trabajos en los que únicamente se seleccionan 
de forma óptima las posiciones de muestreo 
y posteriormente los que incluyen también la 
selección óptima de frecuencias o tiempos de 
muestreo. 


Marco estadístico 


Rouhani (1985) propone el método de 


reducción de la varianza para aumentar 
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una red de monitoreo existente. La red de 
monitoreo se aumenta agregando muestras de 
pozos en lugares que minimizan la varianza 
del error de la estimación de una variable, 
calculada por medio de kriging. Este trabajo 
no involucra la frecuencia temporal de 
muestreo en el diseño. La aplicación de la red 
de monitoreo se hace para la concentración 
de un contaminante, sin embargo, el mismo 
método se puede utilizar en el diseño de redes 
de monitoreo piezométrico. 

Samper y Carrera (1990) presentan una 
comparación de diferentes métodos para 
el diseño de redes de observación basados 
en la implementación del kriging. Dividen 
el diseño de las redes de monitoreo en dos 
grandes categorías: cuando se añade un punto 
de medición a una red ya existente y cuando 
se añaden varios puntos. Incluimos aquí 
únicamente lo referente a la segunda categoría, 
que es la más relevante para este trabajo. En 
este caso se consideran varios métodos para 
añadir posiciones de puntos de medición: 


+ Método de inclusiones sucesivas: se van 
probando uno por uno puntos alternativos 
que se unen a un conjunto dado de pozos y 
se escoge el que minimiza la varianza del 
error de la estimación. 

e Intercambios sucesivos óptimos: se inter- 
cambian los elementos de un conjunto dado 
de puntos existentes con los elementos de 
un conjunto de puntos alternativos para 
encontrar la mínima varianza. 

e Enumeración total: se analizan subcon- 
juntos del total de pozos a monitorear; se 
empieza con todos los subconjuntos de 
un elemento, luego los que se formen con 
dos elementos y así sucesivamente hasta 
llegar al subconjunto que proporcione la 
incertidumbre deseada. 

e ”Rama y límite”: inverso al anterior, el 
conjunto inicial corresponde al conjunto 
total de puntos disponibles y se van su- 
primiendo puntos. 


Dentro de esta categoría también se consi- 
dera la definición de redes de monitoreo 


mediante el uso de mallas regulares que pueden 
tener elementos con diversas geometrías (se 
ubican pozos de una red de monitoreo en cada 
uno de los elementos o nodos de una malla). 
En este caso se trabaja únicamente con los 
parámetros que definen la malla regular y se 
establece una relación entre la varianza y el 
costo asociado con estos parámetros. 

Del análisis realizado para el diseño de una 
red se concluye que el método de inclusiones 
sucesivas, aunque es subóptimo, da buenos 
resultados en la elección de puntos de estima- 
ción, además de que resulta ser una metodolo- 
gía muy económica. Es por lo anterior que esta 
metodología se ha utilizado frecuentemente, 
obteniendo buenos resultados. 

Cameron y Hunter (2000) propusieron un 
método para reducir la redundancia espacial 
y temporal por medio de dos algoritmos: uno 
temporal y otro espacial. El algoritmo temporal 
combina series de tiempo de datos de muchos 
pozos para construir un semivariograma 
temporal compuesto y usa el semivariograma 
para definir frecuencias de muestreo que 
proporcionen datos temporales sin correla- 
ción. En el algoritmo espacial se obtiene una 
estimación inicial de la pluma del contaminante 
usando kriging con todos los pozos existentes. 
Se asignan pesos numéricos, llamados pesos de 
kriging global, a la posición de los pozos, para 
evaluar su contribución relativa a la estimación 
de la pluma. Subconjuntos de pozos con los 
pesos de kriging global más bajos se remueven 
temporalmente y se hace una nueva estima- 
ción de la pluma. La variancia del kriging de 
la nueva estimación se compara con la de la 
estimación inicial para ver si la incertidumbre 
ha aumentado substancialmente. Si no es 
así, se quita el subconjunto de los pozos; en 
caso contrario, el subconjunto de pozos no se 
remueve. En este trabajo no se considera en el 
análisis la correlación cruzada espacio-tiempo 
de la concentración de contaminante. 

Nunes et al. (2004) buscan optimizar 
redes de monitoreo del agua subterránea, 
considerando una reducción en la redundancia 
espacial y/o temporal. Se proponen tres 
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modelos de optimización para seleccionar el 
mejor subconjunto de estaciones de una red 
de monitoreo del agua subterránea: 1) uno 
que maximiza la certidumbre espacial; 2) uno 
que minimiza la redundancia temporal, y 
3) un modelo que maximiza la certidumbre 
espacial y minimiza la redundancia temporal. 
Los modelos de optimización propuestos se 
resuelven con recocido simulado, junto con un 
algoritmo de parametrización, usando entropía 
estadística. Los tres modelos se derivan de una 
ecuación que considera dos términos: uno 
espacial y otro temporal; los modelos emplea- 
dos en el artículo resultan de simplificaciones 
a la ecuación objetivo. La ecuación general de la 
que se derivan los modelos contiene un término 
de varianza y un término que considera series 
de tiempo; las series de tiempo son representa- 
das por funciones matemáticas comunes 
basadas en juicio empírico y experiencia. 

Kumar et al. (2005) diseñaron una red de 
monitoreo de la carga hidráulica basándose 
en la teoría de las variables regionalizadas. 
Se le quita la tendencia a los datos al ajustarle 
un plano. Se ajusta un modelo lineal al 
semivariograma de los residuos obtenidos. Se 
estiman valores mediante kriging universal 
y se determinan los errores en la estimación 
mediante un procedimiento de validación 
cruzada (jacknife). 

Faisal et al. (2007) optimizaron una red 
de monitoreo de la carga hidráulica para la 
estimación de un balance de aguas subterrá- 
neas en una cuenca de la India. Para la 
determinación del balance se utilizó el método 
de la fluctuación del doble espejo de agua 
(llamado así porque considera época de secas y 
de lluvias) propuesto por Maréchal et al. (2006). 
Se llevó a cabo un análisis geoestadístico 
para los datos recabados de carga hidráulica. 
El criterio de optimización de la red se basó 
en la validación cruzada del modelo de 
semivariograma propuesto. Se le asigna un 
orden de importancia a cada pozo de acuerdo 
con el tamaño del error promedio de tres años, 
resultado de la validación cruzada en cada 
sitio. Se realizó el cálculo del balance de aguas 


subterráneas para distintos subconjuntos de 
pozos (siguiendo el orden de importancia 
obtenido de la validación cruzada). Se seleccio- 
nó el número de pozos suficiente para no tener 
pérdida de información en la determinación 
del balance del agua subterránea. 

Para el rediseño de una red de monitoreo de 
los niveles del agua en el acuífero de la Comarca 
Lagunera, México, Esquivel-Victoria (2008) 
utilizó un método geoestadístico propuesto 
por Ahmed (2004) para el diseño óptimo de 
redes de monitoreo de la temperatura del 
aire. En este método, el área de estudio se 
divide en una malla regular, y se obtienen la 
estimación y desviación estándar del error 
de la estimación para cada bloque utilizando 
kriging ordinario para los bloques con el total 
de datos existentes. La optimización se basa en 
la comparación de la desviación estándar del 
error de la estimación sobre los bloques, con un 
valor máximo preestablecido por el usuario, 
conforme se van quitando o añadiendo datos 
en la estimación. 


Marco de modelación 


Al parecer, Loaiciga (1989) fue el primer 
autor que propuso un método para el diseño 
espacio-temporal de redes de monitoreo del 
agua subterránea en el que la redundancia 
espacial y temporal de una red de monitoreo 
se analizaron en conjunto. El método combina 
kriging con una ecuación de transporte 
estocástica para predecir la incertidumbre 
del error de las concentraciones de un 
contaminante estimadas. Él usó el método de 
kriging generalizado para incluir el dominio 
en espacio y tiempo, y de esa manera evaluar 
cómo una muestra tomada de un pozo dado 
en un tiempo dado reduce la incertidumbre 
de todas las posiciones y tiempos en los 
que se obtendrían estimaciones. Para poder 
resolver el sistema de kriging es necesario 
contar con la matriz de covariancia espacio- 
temporal del error de la concentración en las 
posibles posiciones y tiempos de muestreo. 
Loaiciga deriva los elementos de la matriz 


Briseño-Ruiz et al., Método para el diseño óptimo de redes de monitoreo de los niveles del agua subterránea 


de covariancia de la ecuación de advección- 
dispersión que gobierna el transporte de masa, 
discretizada numéricamente. Usó un término 
aditivo de error aleatorio en la ecuación de 
transporte que considera la incertidumbre 
del modelo; este término cambia de un paso 
de tiempo a otro en el modelo, y en cada 
nuevo paso es independiente de todos los 
errores previos. Para la optimización de la red, 
Loaiciga busca determinar dónde y cuándo 
muestrear, para minimizar la variancia del 
error de la concentración estimada (obtenida 
del sistema de kriging espacio-temporal) en las 
posiciones de interés, sujeto a constricciones 
de costo e insesgo. La determinación de un 
plan de muestreo óptimo se plantea como un 
problema de programación entera mixta. 
Yangxiao et al. (1991) combinaron el pro- 
cedimiento de estimación de parámetros 
propuesto en un trabajo previo (Van Geer et 
al., 1991) con un método para el diseño de 
redes de monitoreo de los niveles del agua 
subterránea. El objetivo del diseño de la red 
de monitoreo es observar los cambios espacio- 
temporales de los niveles del agua causados 
por la extracción de agua subterránea. Se 
utiliza un modelo de flujo determinista, al que 
se le suma un error incierto con correlación 
espacial. Los autores proponen utilizar un 
Filtro de Kalman para estimar los parámetros 
de la ecuación de flujo (llamados parámetros 
deterministas) y los del modelo de la matriz 
de covarianza del error (llamados parámetros 
estocásticos). La calibración se realiza para un 
periodo en el que las condiciones estacionarias 
se cumplen y en la que todas las matrices 
necesarias en el algoritmo del Filtro de 
Kalman no cambian con el tiempo. El método 
para la estimación de los parámetros se aplica 
de forma secuencial: primero se propone 
una estimación del modelo de la matriz de 
covarianza del error y la primera ronda de 
calibración de los parámetros deterministas 
se lleva a cabo; estos parámetros se mantienen 
constantes y se obtiene la primera ronda de 
calibración de los parámetros estocásticos. 
Estos dos pasos se repiten varias veces hasta 


que un criterio de error preestablecido para 
ambos conjuntos de parámetros se cumple. 
Se analizaron dos ejemplos de aplicación para 
el diseño de redes de monitoreo; en ambos, 
las frecuencias de muestreo se mantuvieron 
fijas y la densidad de la red se minimizó 
con la restricción de un valor umbral para la 
desviación estándar del error de la estimación. 
Se analizaron varias alternativas y la mejor fue 
elegida por inspección. Los autores discuten 
la importancia relativa de la densidad de 
la red espacial y la frecuencia de muestreo 
en relación con el tiempo de respuesta 
del sistema. Encontraron que si el sistema 
reacciona rápido, la optimización espacial de 
la red es importante. Si el sistema reacciona 
lentamente, tanto la optimización temporal 
como la espacial son importantes. 

Herrera (1998) propuso una metodología 
para el diseño óptimo espacio-temporal de 
redes de monitoreo de la calidad del agua 
subterránea. El método combina un Filtro de 
Kalman estático con un modelo de transporte 
estocástico para predecir la varianza del error 
Una de 
las entradas para el Filtro de Kalman es la 


de las concentraciones estimadas. 


matriz de covarianza espacio-temporal de la 
concentración del contaminante en todas las 
posiciones y tiempo de muestreo posibles. 
Los elementos de esta matriz de covarianza se 
derivan de la ecuación de un modelo numérico 
de transporte de masa con advección y dis- 
persión que gobierna el transporte de masa 
numéricamente discretizada. En este modelo, 
las variables aleatorias que se consideran son 
la conductividad hidráulica y la fuente del 
contaminante. Para optimizar la red se busca 
determinar dónde y cuándo muestrear para 
minimizar una función de la varianza del 
error de la concentración estimada (obtenida 
del Filtro de Kalman) en las posiciones y 
tiempos de interés. La determinación de un 
plan de muestreo óptimo se obtiene a través de 
un procedimiento secuencial, que selecciona 
los puntos de muestreo en espacio-tiempo 
que minimizan la función en cada paso, y se 
detiene cuando la función alcanza un valor 
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predeterminado. Este método de optimización 
se puede considerar como uno de inclusiones 
sucesivas, que incluye tanto la selección de 
las posiciones de los pozos como los tiempos 
de monitoreo. En este trabajo se presentaron 
algunos ejemplos sintéticos de aplicación. 

Posteriormente, Herrera et al. (2001) aplica- 
ron la metodología propuesta por Herrera 
(1998) en el diseño de una red de monitoreo 
de la calidad del agua para el acuífero Toms 
River en New Jersey, Estados Unidos. Y en 
Herrera y Pinder (2005) incluyeron, además 
de la varianza, al coeficiente de variación 
como una medida de la incertidumbre. Como 
se mencionó antes, la metodología utilizada 
en el presente trabajo es una modificación de 
la propuesta por Herrera y Pinder. 

Zhang et al. (2005) diseñaron una red de 
monitoreo óptima de la calidad del agua sub- 
terránea, combinando un Filtro de Kalman y un 
algoritmo genético para reducir al máximo el 
coeficiente de variación en posiciones y tiempos 
establecidos. La metodología propuesta es una 
modificación de la de Herrera y Pinder. En este 
caso, para diseñar redes de monitoreo de la 
calidad del agua: a) se generan realizaciones de 
la conductividad hidráulica usando la técnica 
Latyn Hypercube Sampling (LHS), suponiendo 
a la conductividad como un proceso aleatorio; 
b) se simula el flujo y transporte del agua 
subterránea para generar realizaciones de la 
concentración del contaminante, usando las 
realizaciones de la conductividad obtenidas 
en el paso anterior; c) se calcula una matriz 
de covarianza espacio-temporal de las con- 
centraciones del contaminante; d) se usa un 
algoritmo genético para seleccionar posiciones 
y tiempos de muestreo que reduzcan el coe- 
ficiente de variación al máximo, en posiciones 
y tiempos de estimación. Se usa un Filtro 
de Kalman para actualizar la matriz de 
covarianza espacio-temporal después de que 
se ha tomado una muestra; e) se sigue con el 
proceso de selección hasta que se logran los 
valores objetivo del coeficiente de variación y 
se alcanza un escenario óptimo. Se realizó un 
problema ejemplo en Tomas River, New Jersey, 
Estados Unidos. 


Hasta donde sabemos, no existen trabajos 
publicados en los que, al igual que en este tra- 
bajo, se utilicen modelos de flujo subterráneo 
que incluyan la incertidumbre de parámetros 
con correlación espacial para diseñar redes 
de monitoreo óptimas de los niveles del agua 
subterránea, ya que Herrera (1998) mostró que 
el modelo que utilizan Yangxiao et al. (1991) 
no es adecuado para representar los errores 
obtenidos como efecto de estos parámetros 
con correlación espacial. 


Metodología de diseño óptimo espacio- 
temporal de una red de monitoreo 
piezométrica 


El método propuesto originalmente por He- 
rrera (1998) para el diseño óptimo de redes 
de monitoreo de calidad del agua subterránea 
se adaptó y probó en este trabajo para la 
optimización de redes de monitoreo de los 
niveles del agua subterránea. El método 
que se propone es adecuado para optimizar 
redes de monitoreo de largo plazo, en casos 
en los que el acuífero ha sido investigado 
ampliamente y se cuenta con un modelo de 
flujo para el mismo. La metodología se divide 
en dos procedimientos: 1) uno que se usa para 
estimar la carga hidráulica (en lo siguiente )) 
y la incertidumbre del error de esta estimación 
cuando se tienen datos de h en diferentes pozos 
y tiempos de monitoreo, y 2) un método que 
escoge las posiciones de los pozos y los tiem- 
pos de muestreo que minimicen la predicción 
de la incertidumbre de h obtenida en el paso 1, 
con los cuales se definen la red de monitoreo 
y su programa de muestreo. En los párrafos 
siguientes se describen con mayor detalle las 
dos etapas del método. 


1. Estimación de la carga hidráulica y de la 
incertidumbre de su error. Esta parte del 
método emplea un Filtro de Kalman y un 
modelo estocástico de flujo para obtener 
una estimación de h y la incertidumbre 
que el error de ésta tendría si datos de h 
de muestras tomadas en un conjunto de 
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posiciones en diferentes tiempos se usaran 
para estimar la carga hidráulica. Primero se 
emplea el modelo estocástico de flujo para 
calcular una estimación inicial (o a priori) 
de h en diferentes posiciones y tiempos, 
así como su matriz de covarianza espacio- 
temporal. El modelo estocástico tiene la 
ventaja de que sintetiza el conocimiento 
que se tiene del acuífero y a la vez permite 
reconocer áreas con diferentes niveles de 
incertidumbre (a través de la varianza). 
Posteriormente, con el Filtro de Kalman 
se actualizan tanto la estimación inicial 
como la matriz de covarianza, al incluir 
los datos de piezometría disponibles y las 
posiciones de los pozos correspondientes. 
Es importante recalcar que para obtener la 
incertidumbre del error de la estimación no 
es necesario conocer el valor de los datos, 
sino únicamente las posiciones de las que 
se obtendrían los mismos; en este sentido, 
es posible predecir esta incertidumbre sin 
haber tomado las muestras. 

2. Selección óptima de las posiciones de los 
pozos y los tiempos de muestreo. Una 
función de la varianza del error de la 
estimación, calculada con el procedimiento 
1, se utiliza como criterio para escoger 
los pozos de muestreo y los tiempos de 
muestreo en cada pozo. La función que se 
utiliza depende de los objetivos del diseño. 
Para minimizar la función se pueden utilizar 
diferentes métodos de optimización; hasta 
ahora se ha usado un método secuencial 
que escoge los puntos espacio-temporales 
que la minimizan en cada paso y se detiene 
cuando el usuario lo decide. 


Modelo estocástico 


La ecuación de flujo se emplea para describir la 
evolución de h: 


Y (K-Vh)=55 -+Q=0 (1) 


en donde h es la carga hidráulica [L]; K es la 
conductividad hidráulica [L/T]; Q representa 


fuentes o sumideros [L/T], y Ss es el coefi- 
ciente de almacenamiento específico [1/L]. 
Esta ecuación se  discretiza usando un 
esquema numérico y el sistema se resuelve 
para el vector h de las incógnitas de la carga 
hidráulica en un subconjunto de puntos 
asociados con la malla numérica. A esta malla 
la llamaremos malla del modelo estocástico. 
Cualquiera de los parámetros, así como las 
condiciones iniciales o las de frontera, pueden 
representarse como variables aleatorias. En 
general se reconoce que uno de los parámetros 
con mayor incertidumbre en esta ecuación es 
la conductividad hidráulica. Por tal razón, en 
el ejemplo presentado en este trabajo, la con- 
ductividad hidráulica se representa como un 
campo aleatorio espacialmente correlaciona- 
do. Para caracterizar estadísticamente a K 
se calcula su media y se realiza un análisis 
geoestadístico, a fin de obtener un modelo de 
semivariograma que representa su estructura 
espacial. 


Filtro de Kalman 


El Filtro de Kalman es un conjunto de ecuaciones 
matemáticas que obtienen recursivamente 
estimaciones lineales no sesgadas y con va- 
rianza mínima para el estado de un sistema 
utilizando datos con ruido (Jazwinski, 1970). 
Asimismo, se puede emplear para actualizar la 
estimación en el momento en que se cuente con 
datos nuevos de la variable, sin necesidad de 
utilizar datos anteriores. El término recursivo 
significa que el filtro recalcula la solución cada 
vez que una nueva observación o medida es 
incorporada en el sistema. 

En su forma general, el filtro se basa 
en dos ecuaciones: una ecuación dinámica 
y una ecuación de medición. La ecuación 
dinámica contiene un error aleatorio que no 
tiene correlación en el tiempo. Sin embargo, 
Herrera (1998) mostró que cuando se usa un 
campo aleatorio con correlación espacial para 
representar a la velocidad y se deriva con 
base en éste el modelo de transporte como 
una ecuación con coeficientes deterministas 
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más un error aleatorio, el error puede tener 
una fuerte correlación en el tiempo. Por esta 
razón, el método que propuso Herrera (1998) 
emplea lo que llamamos el Filtro de Kalman 
estático, que únicamente emplea la ecuación de 
medición, incorporando el tiempo a través de 
la utilización de vectores espacio-temporales. 

La ecuación lineal de mediciones del Filtro 
de Kalman discreto, que relaciona al vector h 
de la variable en las posiciones y tiempos en 
que se quiere estimar con las muestras z es: 

2¡=Hj¡h+0; (2) 

en donde tz, j=1,2,...) es una secuencia de 
mediciones de la carga hidráulica. La matriz 
de A, es una matriz de 1 * N, que no es cero 
únicamente en la posición correspondiente a la 
entrada de h, en donde se toma el muestreo j, y N 
es la dimensión del vector h. H,es la matriz del 
muestreo j; h = (h,,) es el vector espacio-tiempo 
con los valores estimados de la carga hidráulica 
en las posiciones y tiempos de interés (5, es la 
carga hidráulica en la posición x,en el tiempo 
£,); (o, j=1,2,...) es un escalar que representa 
el error de medición. Esta es una secuencia 
blanca Gaussiana, con media cero y varianza 
1, La secuencia del error de medición to) y el 
vector h son independientes. 

La matriz de covarianza del error de la 
estimación es: 


PE ins) 


en donde hr" = El (h/z,, Zy...,2,)"] es el valor 
esperado de h, dadas las mediciones z,, 


(3) 


Dis pZ 
notación, el superíndice identifica el número 


y T denota transpuesta. En esta 


n de mediciones que se usan para obtener la 
estimación. 

En la aplicación del Filtro de Kalman, 
se estima h en una submalla de la malla del 
modelo estocástico, a la que llamamos malla 
del Filtro de Kalman. Esta malla tiene dos tipos 
de nodos, definidos como las posiciones de los 
posibles pozos de muestreo (S) y de los puntos 
de estimación (E), que se definirán cuando se 


explique el método de optimización de la red 
de monitoreo. Para la implementación del filtro 
se requiere proponer una estimación a priori 
de h en espacio y tiempo (h1%), así como de la 
matriz de covarianza del error de estimación 
(P%). Dadas estas estimaciones a priori, la 
estimación lineal de mínima varianza para h 
se puede obtener secuencialmente a través de 
las siguientes fórmulas: 


pa = p" + Koa (2 0 Ars ñ”) (4) 


Pri = PO Kris Has P" (5) 


1 
Ka =P "Hal HnaP"Hia tra) (6) 


Estimación de los momentos a priori 


El procedimiento para calcular la estimación a 
priori de hen espacio y tiempo (h0), así como 
la matriz de covarianza del error de esta 
estimación (P%) es la simulación estocástica 
o Monte Carlo (Zhang, 2002), y se explica a 
continuación. 


1) Dadas la media y la de covarianza (o 
semivariograma) de la conductividad 
hidráulica se obtienen N, realizaciones 
de esta variable que sean acordes con las 
mismas. 

2) Se resuelve numéricamente la ecuación 
de flujo empleando cada realización, y se 
obtienen realizaciones de h para un periodo 
de tiempo. 

3) Utilizando promedios de las realizaciones 
de h para cada nodo de la malla del 
Filtro de Kalman, obtenemos la media en 
espacio-tiempo, la que se propone como 
la estimación a priori, esta estimación es 
un vector al que denotaremos /* y de la 
misma manera se calcula su matriz de 
covarianza. La estimación a priori de h en 
espacio y tiempo (/1), así como la matriz de 
covarianza del error de la estimación (P0), 
se obtienen con las siguientes ecuaciones: 
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1 Nr 


hy == Y hi (7) 
p N, = p 


Cov(e;, ig )> Covip:jg 


1 Nr 


= N, =1 2 Hp Js, ha) (8) 


en donde hi denota la k-ésima realización 
de h en la posición x,en el tiempo b h" es el 
vector que contiene la media h,, para todas las 
realizaciones he; Cove, y e,,) es el elemento (i, p; 
j, 9) de la matriz de covarianza del error de la 
estimación P%, donde el error e = hi, - h.; como 
h" es un vector en espacio-tiempo, éste incluye 
una estimación de h para todas las posiciones y 
tiempos de interés; por tanto, Pes también una 
matriz de covarianza en espacio-tiempo. 

Para el ejemplo que se presentará, el núme- 
ro total de entradas que integran la matriz 
de covarianza del error de la estimación de 
la carga hidráulica es Nt = (m*p)?, donde p es 
el número total de tiempos de estimación y 
posibles tiempos de muestreo, y m es el número 
de nodos de la malla del Filtro de Kalman; m 
= NE + NS, donde NE es el número de puntos 
de estimación y NS es el número de pozos de 
muestreo, ya que no existe intersección entre 
los puntos de estimación y los de muestreo. 


Función para minimizar la varianza del 
error de la estimación 


Como se explicó antes, para la selección Óptima 
de las posiciones de los pozos y los tiempos de 
muestreo se utiliza un método de optimiza- 
ción, en el que se minimiza una función de 
la varianza del error de estimación. En este 
caso, esta función es la suma de la varianza 
del error sobre todas las posiciones y tiempos 
de estimación. A esta función la llamamos 
la varianza total de la estimación de la carga 
hidráulica y se denota por 0?,(11). Ésta, cuando 
se calcula con el plan de monitoreo que incluye 
n muestras, es: 


oñ(n) = Y oj(n) (9) 


i,peE 


en donde 0%, (11) es la varianza del error de la 
estimación en el ¡-ésimo punto de estimación y 
en el tiempo de estimación p, y E es el conjunto 
de los índices de los puntos de estimación. 
Las varianzas 0%, (1) se obtienen del Filtro de 
Kalman después de tomar n muestras, y son 
los elementos de la diagonal de la matriz de 
covarianza P", que corresponden a los puntos 
y tiempos de estimación. 

El algoritmo iterativo de optimización 
para la selección de las posiciones y tiempos 
de muestreo se explica diagramáticamente 
en la figura 1. La selección de las posiciones 
y tiempos de muestreo se hace en forma 
secuencial. El método utiliza un algoritmo 
iterativo de optimización que escoge una a 
la vez las posiciones de los pozos y las fechas 
de muestreo que minimizan la varianza total 
de la estimación, calculada como la suma de 
la varianza sobre todos los nodos y tiempos 
en que se va a estimar la h. En otras palabras, 
dadas las posiciones x,,xX,,...,Y, de los posibles 
pozos de muestreo, y la matriz de covarianza a 
priori PY, en el primer paso del procedimiento, 
para x, 1 = 1,...s, con el Filtro de Kalman se 
calcula la varianza del error de la estimación 
que se obtendría al agregar el dato de cada 
posición de muestreo y se escoge la posición 
y; que da la varianza total menor. A la matriz 
de covarianza resultante de aplicar el filtro 
con P? y yí la denotamos P;. Una vez que ya se 
seleccionaron las posiciones de los pozos, Yi 
j=1...,k para formar parte de la red óptima, se 
toman las posiciones que no han sido escogidas, 
y con la matriz de covarianza actualizada P;, 
que incluye los cambios debidos al muestreo de 
las Yi ¡ = 1,...,k posiciones anteriores, se aplica 
el Filtro de Kalman y nuevamente se escoge 
la posición que da la varianza total menor. 
Las fórmulas empleadas para minimizar la 
varianza total se presentan en Herrera (1998), y 
Herrera y Pinder (2005). 

Debido a que el proceso es secuencial, se 
realiza un análisis de la contribución en la 


Tecnología y 


Ciencias del Agua, vol. 1H, núm. 4, octubre-diciembre de 2011 


Briseño-Ruiz et al., Método para el diseño óptimo de redes de monitoreo de los niveles del agua subterránea 


Con método llamado 
simulación secuencial 


Gaussiana de librerías 


Análisis geoestadístico 
GSLIB, generar las 


Calcular la media del LnK y el semivariograma 
realizaciones de LnK 


ES Obtener realizaciones de LnK | 


Generar modelo estocástico 


Resolver la ecuación de flujo 


v(k-v1)- ss E +0=0 


Con los promedios de las cargas hidráulicas en cada posición y tiempo, 
obtener la media y su matriz de covarianza 


Le 
COV jp: jq N 20 1p 
r k= 


Dada la estimación a priori ho = (1) y 


la matriz de covarianza a priori P0= (Cow, in) 


Aplicar el filtro de KALMAN para: 


2. Aplicar el filtro y calcular las posiciones 
subsecuentes de la red 


y 


Con la matriz de covarianza P*! como la 
matriz a priori PO, calcular las varianzas 
al agregar las posiciones X,.+2..hasta s 


P?=P%* calcular PY =P" -K ¿HP 


hasta calcular las matrices de covarianza de 
las posiciones restantes de muestreo 


A 
I 


1. Seleccionar la posición del pozo X;,,1 
de los pozos de muestreo 


Calcular la varianza de la estimación de la 


carga hidráulica que se obtendría al agregar 


la posición x,,,1 con hy y PO sl sl 


prtl =p0 -—K.H po 


n+12n41 
en donde 


-1 
0pyT 07yT 
Kara =P Hsa(Hps1P Ha + 2) 


Método secuencial 


Seleccionar una posición, que minimice la 
función en cada paso 


De 2 1 s1? 
o7(s1+1) =07(s1) pe Y (25) 


Método secuencial 


hasta calcular las matrices de covarianza 


(hasta s pozos de muestreo) 


para cada posición de muestreo 


Seleccionar la posición que posea la 
varianza total menor, 


Establecer la matriz de covarianza P*!*" 


seleccionada, como la matriz a priori, y 
calcular las varianzas para las posiciones 


on) = Yo; (1) 


¡,peE 


Xn+3... hasta las s posiciones totales... 


y su matriz de covarianza P* (espacio- 
tiempo). La posición seleccionada es la 


Detener cuando se alcanza el valor 


primera que conforma la red de muestreo pretederminado de la varianza del error 


Postprocesamiento 
po 
ho,P y ltvtata 


Con el Filtro de Kalman actualizar la 
estimación 


Ln+l a Gn 
h =h" + Rot na _ H, +12 ) 


K =P” yl H,,,p” yl E 
a mi tn+1 14 Th 41 


Figura 1. Método de Herrera y Pinder. 


a] 


(hasta s pozos de muestreo) 


iencias del Agua, vol. Il, núm. 4, octubre-diciembre de 2011 


E 


« vJ30]0u32], 


Briseño-Ruiz et al., Método para el diseño óptimo de redes de monitoreo de los niveles del agua subterránea 


reducción de la varianza de cada pozo para 
determinar cuándo parar el proceso y así quede 
definido el número total de pozos de la red de 
monitoreo óptima. 


Caso de estudio sintético 


La metodología descrita se probó en un caso 
sintético basado en una representación sim- 
plificada del acuífero del Valle de Querétaro 
(AVQ). Como una primera prueba del método, 
se rediseñó una red de monitoreo de los 
niveles del agua subterránea, suponiendo que 
los cambios en los mismos son insignificantes. 
Para esto se utilizó un modelo numérico de 
flujo en estado estacionario. 

El AVO está localizado en la porción sur 
oriental del estado de Querétaro, México. Se 
ubica en el subsuelo de la ciudad de Querétaro 
y su zona conurbada (ZMCQ) (Mendoza y 
Herrera, 2007). Su población es la mayor de 
todo el estado, estimada en 962 240 habitantes 
en 2007 (Urquiza, 2008). La mayor explotación 
de agua subterránea en el AVQO es para el uso 
público urbano, con la cual se abastece la 
ZMCO. 


Objetivo de la red de monitoreo 


El objetivo del diseño de la red de monitoreo 
de los niveles del agua subterránea consistió 
en seleccionar, de entre 33 pozos preexistentes 
en el AVO, aquellos que no dieran información 
redundante para la estimación de h en todo 
el acuífero. De acuerdo con este objetivo se 
estimará h en todo el acuífero, dándole el 
mismo peso a todas las zonas. A continuación 
se explica cómo se aplicó el método. 


Simulación secuencial Gaussiana 
Con base en un análisis geoestadístico de 46 
datos de K existentes en la zona de estudio, 


se ajustó un modelo exponencial para el 
semivariograma muestral de Ln K: 


3d 
Yink(d)=c9 +01 (1 o) (10) 


en donde y LnK es el semivariograma de LnK; 
Cy la pepita o nugget; c,+c, el sill o la meseta del 
semivariograma, d, la distancia, y a es la escala 
de correlación. Los valores ajustados para los 
parámetros fueron C,= 0.43, c,= 0.57, a= 8 900 
m y la distancia máxima considerada es de 
1 262.05 metros. 

Con el modelo ajustado mediante el 
análisis geoestadístico se generaron cuatro mil 
realizaciones aleatorias de Ln K mediante el 
método conocido como simulación secuencial 
Gaussiana (SGSim) de las librerías del GSLIB 
(Deutsch y Journel, 1998). Para generar las 
simulaciones mediante SGSim se requiere quec, 
+c,= 1;en el modelo ajustado esto se cumple. La 
malla sobre la que se obtienen las realizaciones 
de LnK tiene elementos cuadrados y cubre en 
su totalidad el área del modelo determinista. 
Mediante una transformación de los valores 
de estas realizaciones se obtienen realizaciones 
de K con media 1.63 m/d (los detalles de 
la transformación se pueden consultar en 
Manchuk et al., 2004), que es igual a la media 
de los datos de K disponibles en este acuífero. 


Modelo determinista y estocástico de 
flujo de agua subterránea 


La simulación de flujo se realizó con un simula- 
dor determinista convencional de elemento 
finito llamado Princeton Transport Code (Babu 
et al., 1993). Éste se adaptó utilizando el modelo 
calibrado por Simuta (2005). Las características 
del modelo del AVQ se pueden ver en la 
figura 2. La zona de modelación corresponde 
a la del acuífero del Valle de Querétaro e 
incluye también al Valle de Obrajuelo, en el 
estado de Guanajuato, con un total de 426.87 
km?. Este último acuífero se incluye debido 
a que aunque por razones administrativas 
la Comisión Nacional del Agua incluye al 
Valle de Obrajuelo como parte del acuífero 
de Celaya, Simuta (2005) concluye que existe 
conexión hidráulica entre ellos y es necesario 
modelarlos en forma conjunta. La malla del 
modelo es triangular y cuenta con 6 676 nodos, 
13 066 elementos y una capa. Los valores de 
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conductividad hidráulica en cada realización 
se asignaron a los nodos de la malla numérica 
del modelo determinista, encontrando el nodo 
más cercano al centro de la malla de elemen- 
tos cuadrados utilizada para las realizaciones 
de K. Los parámetros hidráulicos del modelo 
son: el almacenamiento específico es uniforme 
en toda la zona de modelación (Ss = 0.001). La 
recarga se considera variable en la zona con 
rangos de 0.011 a 0.067 metros/día para la 
recarga vertical originada por la precipitación y 
de 2.45 e* a 0.0003229 metros / día en la recarga 
originada por el riego. En el caso de estudio no 
se considera extracción por bombeo. Se tienen 
cinco zonas con condiciones de frontera de 
carga asignada, en promedio cien metros por 
debajo de la topografía. 

El modelo estocástico de flujo tiene la misma 
malla que el modelo determinista de flujo. El 
diseño de la malla cumple el criterio sugerido 
por Ababou (1998) para que las simulaciones 
del LnK no tengan una distorsión significativa 
debido a la discretización. Éste es, Ax< 4 
donde Ax es el tamaño de los elementos de L 
malla y A es la escala de correlación que para 
el modelo exponencial es 2/3. Las condiciones 
de frontera, así como todos los parámetros 
del modelo se tomaron como deterministas, 


con excepción de K, que es considerada 
como un campo aleatorio con media de 1.63 
metros /día y correlación espacial dada por 
el semivariograma exponencial descrito en 
la sección anterior. Para obtener la solución 
del modelo estocástico se generaron cuatro 
mil realizaciones de K y sus correspondientes 
realizaciones de h. 

Las posiciones en las que se estimará h 
con los datos que se obtendrán de la red de 
monitoreo se asocian con los nodos de lo que 
llamamos malla de estimación. De acuerdo 
con lo explicado anteriormente, en estos 
puntos de estimación es en donde se minimiza 
la varianza del error de estimación de h en el 
proceso de optimizar la red de monitoreo. De 
acuerdo con el objetivo de la red de monitoreo, 
para darle el mismo peso a toda el área del 
acuífero en la optimización se propuso como 
malla de estimación una malla compuesta por 
66 puntos de estimación aproximadamente 
equi-espaciados sobre el AVO (figura 3). Las 
estimaciones de h se calcularon en los nodos de 
la malla de estimación, y por lo tanto la matriz 
de covarianza incluyó únicamente estos nodos 
y las posiciones de los 33 posibles pozos de 
monitoreo, que también se pueden ver en la 
figura 3. 


Vista en isométrico 


O Valle de Querétaro 


O Flujo en estado 
estacionario 


O No se considera bombeo 


y Almacenamiento específico 
Ss = 0.001 


o Conuctividad hidráulica 
K=1.63 m/día 


O Condiciones de frontera con 


cargas asignada —————8' 


O Recarga por lluvia e 
irrigación 


Malla triangular con 
6 676 nodos 


ó Norte 


Vista en planta 


Figura 2. Modelo determinista. 
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Para verificar que el número de las Resultados 
realizaciones del modelo estocástico fueran 
suficientes, se estableció el criterio de La figura 5 muestra la varianza del error de 
convergencia para la matriz de covarianza la estimación inicial, ésta es la varianza de h 
basado en las siguientes ecuaciones: obtenida de la simulación estocástica. Como 


se puede observar, tiene valores entre 0 y 9 m?, 

yA con valores mayores hacia el límite oeste de la 

Nina) vi yl dif, maayli, j)| (11) zona en la que se obtendrán las estimaciones 
á y en la región noreste también se tiene 

un máximo local. Esto indica que la carga 


difonm+a) (1,5) = Bn (i, J) — Peri) (ls J) (la) 


en donde P, es la matriz de covarianza calcula- 


da con m realizaciones, A es el incremento en el 


número de las realizaciones y N es el número 0.200 


de entradas de la matriz de covarianza. En la 


prueba, el incremento A utilizado fue de 20. Con e 
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mil realizaciones son suficientes (ver figura 4) 
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Figura 3. Pozos de monitoreo (izquierda) y malla de estimación (derecha). 
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Figura 5. Varianza del error de la estimación inicial, y varianza final después de muestrear 12, 22 y 33 pozos. 


hidráulica es más sensible a los cambios en la 
conductividad hidráulica en esas regiones. 

Como se mencionó anteriormente, el 
método de optimización propuesto es se- 
cuencial y se añade un pozo de monitoreo en 
cada paso. Esto produce los resultados que 
se muestran en forma de gráfica en la figura 
6, en donde se representa la varianza total en 
los puntos de estimación contra el número de 
pozos incluidos en la red. Se observa que los 
primeros pozos elegidos proporcionan mucha 
información y la varianza total se reduce 


significativamente; conformese vanañadiendo 
pozos a la red de monitoreo, la información 
obtenida por cada pozo va disminuyendo. La 
varianza total inicial (m?) para los 66 puntos 
de estimación es de cerca de 120, y conforme 
se van añadiendo pozos a la red de monitoreo, 
esta varianza disminuye hasta un valor 
cercano a 40. Esto indica que al utilizar todos 
los pozos de la red de monitoreo, la varianza 
total no baja de este valor, que al promediar 
sobre los 66 puntos de estimación da una 
varianza cercana a 0.6. Sin embargo, esta 
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Figura 6. Varianza total versus número de pozos en la red. 


figura también hace ver que los últimos pozos 
contribuyen muy poco a reducir la varianza, 
por lo que es natural analizar la pertinencia de 
dejarlos en la red de monitoreo. 

Por la forma en la que se construye la red de 
monitoreo, el orden en el que se eligen los pozos 
indica el grado de importancia de cada pozo en 
la disminución de la varianza; por ese motivo, 
a éste le llamamos orden de importancia de los 
pozos. El orden de importancia de los pozos 
se puede ver en la figura 7. Como se puede 
observar, el pozo con orden de importancia 
1 se encuentra en la parte oeste de la zona 
establecida para optimizar la red de monitoreo; 
ello es lógico, ya que ésta es la región con mayor 
varianza. De igual forma, el pozo con orden de 
importancia 2 se ubica al noroeste de la zona de 
estudio, en donde se señaló anteriormente que 
existe un máximo local. Es importante recal- 
car que al seleccionar estos pozos se toma en 
cuenta la varianza del error de la estimación 
en toda la zona de influencia del pozo y no 
únicamente en el punto que se escoge. De esta 
manera, se escogen también los pozos con orden 
de importancia 3, 4 y 5, en posiciones en las que 
se tienen varianzas altas, y en las que no se han 
escogido pozos anteriormente. A partir del 
pozo con orden de importancia 6, van cobran- 
do más importancia las selecciones realizadas 
anteriormente y las zonas de influencia de los 
pozos (zonas con correlaciones altas obtenidas 
de las matrices de covarianza resultantes de 
aplicar el Filtro de Kalman conforme se van 


añadiendo los pozos), y se inicia la selección 
de pozos cercanos a otros ya elegidos o en la 
zona central, que inicialmente tiene varianzas 
pequeñas. 


Análisis para determinar el número total 
de pozos en la red de monitoreo 


Como se comentó antes, en la optimización de 
esta red de monitoreo el objetivo es determi- 
nar los pozos que no sean redundantes. Para 
hacer esto se lleva a cabo un análisis en el 
que se produce una estimación inicial /* y se 
actualiza con los datos asociados con una red 
de monitoreo con un determinado número 
de pozos, mediante la aplicación del Filtro 
de Kalman. En el caso de estudio sintético, 
este procedimiento se realizó con el objeto 
de mostrar cómo se haría este análisis en un 
caso real, y para ello se tomó al azar una de 
las realizaciones de h del modelo estocástico, 
como si fuera la carga hidráulica en campo que 
se quiere estimar. 

Como un primer paso, se muestran los 
valores de h en los 6 676 nodos de la malla 
del modelo (figura 8a). También se muestran 
las estimaciones de h producidas con el Filtro 
de Kalman en los 66 nodos de la malla de 
estimación utilizando 12 (figura 8b), 22 (figura 
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Figura 7. Orden de importancia de los pozos. 
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8c) y 33 (figura 8d) datos de h que correspon- 
den a los pozos con el orden de importancia del 
1 hasta el número indicado. La estimación se 
realizó únicamente para la zona del estado de 
Querétaro (el Valle de Obrajuelo no se incluyó 
en la estimación). La comparación entre estas 
piezometrías indica, por un lado, que la 
estimación generada con el filtro, empleando 
12, 22 y 33 datos de los pozos de muestreo de 
la red propuesta, generan piezometrías acordes 
con la realización escogida de h; por otro lado, 
también las piezometrías generadas con 12, 22 
y 33 datos son muy similares. 


En segundo lugar, se realizó una compara- 
ción entre la varianza inicial del modelo (figura 
5a) y las varianzas finales obtenidas con el 
Filtro de Kalman después de incorporar los 
datos de 12, 22 y 33 pozos de muestreo de la red 
de monitoreo propuesta. Las figuras 5b, 5c y 5d 
muestran las varianzas que se obtienen del Fil- 
tro de Kalman después de muestrear 12, 22 y 33 
pozos. La varianza total después de muestrear 
33 pozos, disminuyó considerablemente y tiene 
valores de entre O y 3. También se puede ver que 
las varianzas del error de la estimación usando 
22 pozos y usando 33 pozos son muy similares, 
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Figura 8. Realización de h y estimación de h con 12, 22 y 33 pozos. 
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lo que confirma que existe una redundancia 
importante entre los 11 pozos que se añaden. 
Adicionalmente, se calcularon los errores 
en los 66 puntos de estimación, restando a la 
realización escogida la estimación producida 
con el filtro. Los errores más grandes detecta- 
dos se encuentran ubicados en las zonas NE y 
SW y varían de -1.29 a 1.93 metros, y se ubican 
en las zonas donde se encuentran las varianzas 
más grandes; ello se debe a que en estas zonas 
no se cuenta con ningún pozo de muestreo 
que aporte datos para reducirlas. Esto indica 
que se requerirían pozos adicionales en tales 
zonas para obtener una mejor estimación. 
Para comparar las alternativas entre la 
opción de una red de monitoreo con 12, 22 o 
33 pozos, se realizó el cálculo de los errores. El 
error e, se calculó mediante la diferencia de los 
datos de la realización de h y de la h estimada 
mediante el Filtro de Kalman en los puntos 
de la malla de estimación con cada una de las 
opciones de números de pozos. Los errores 
calculados son el error medio (EM), el error 
cuadrático medio (ECM) y el error cuadrático 
medio estándar. El EM es un indicador 
de la sobrestimación o subestimación de 
los valores de la realización de h. El error 
cuadrático medio (ECM) es una medida de la 
magnitud del error de estimación. La teoría 
indica que mientras más cercanos estén el 
EM y el ECM a cero, más adecuado será el 
ajuste entre la realización de h y la estimación 
de Í con 12, 22 o 33 pozos (Armando, 1994). 
El error cuadrático medio estándar (ECMEB), 
mientras más cercano esté a 1 indica el nivel 
de confianza entre la realización de h y los 


valores de la estimación h. Los resultados de 
estos errores para cada uno de los casos de 
estudio se muestran en cuadro 1. Se puede 
observar cómo el ECM y ECME tienen 
valores cercanos a lo que indica la teoría, sin 
embargo, es importante notar que cuando 
el número de pozos que conforman la red 
disminuye, también decrece el ECME. 


Red de monitoreo final 


Del análisis anterior se concluye que al 
utilizar una red de 22 pozos se obtiene 
una estimación con errores similares a los 
obtenidos con los 33 pozos, esto indica que 
los once pozos con orden de selección de 23 a 
33 se pueden considerar como redundantes. 
Además, los resultados de la piezometría 
indican que con 22 pozos obtenemos una 
representación piezométrica muy parecida a 
la real, con errores de estimación del orden de 
1.92 a -1.33 metros, y con un error cuadrático 
medio estándar de 0.97, muy cercano a 1. Es 
importante tener en cuenta que el número de 
pozos que conformarán la red de monitoreo 
final estará definida por la autoridad federal, 
estatal o municipal de aguas, o en su caso 
por los contaminadores potenciales del agua 
subterránea, como empresas que requieran 
de algún monitoreo de calidad del agua en 
sus instalaciones. La red está en función del 
tipo y volumen de datos requeridos, y varía 
considerablemente en cuanto al aspecto de 
gestión de que se trate, pero inevitablemente 
también dependerá delos recursos financieros 
disponibles y del nivel de riesgo que se esté 
dispuesto a correr. 


Cuadro 1. Comparación de los errores en la estimación con las opciones consideradas. 


Número de és E de ea E ey Máx error | Máx error 
pozos (en = - 2 e ) [Ecos == 2 (e;) ) ECME= Y 7 A Ps 
33 0.1811047 0.38521572 0.96964912 1.935 1.29 
22 0.18506818 0.406800151 0.927066718 1.943 1.325 
12 0,12124242 0.42495848 0.72998227 2.287 1.46 
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Del análisis de la varianza de los errores 
de las estimaciones, al utilizar todos los pozos 
posibles, indica que en las zonas NE y SW de 
la zona de estudio sería recomendable incluir 
algunos pozos de monitoreo que aporten datos 
de carga hidráulica y con ello se logre una 
mejor estimación de h en la zona de estudio. 


Conclusiones 


Del análisis anterior se concluye que el diseño 
de la red de monitoreo piezométrica para este 
caso de estudio mediante el método de Herrera 
y Pinder es satisfactorio, y que se tienen buenas 
perspectivas para aplicar exitosamente esta 
metodología al diseño de redes de monitoreo 
de los niveles del agua subterránea en casos 
reales. 

Como se mostró en los resultados del aná- 
lisis de la varianza en la figura 6, el mode- 
lo estocástico de flujo es una valiosa he- 
rramienta para evaluar la incertidumbre en 
las estimaciones del modelo, haciendo la 
metodología muy efectiva en la selección de 
las posiciones y los tiempos de muestreo que 
minimizan la incertidumbre de la estimación. 
Asimismo, la metodología es capaz de sugerir 
en qué zonas se requiere la construcción y/o 
habilitación de pozos de monitoreo que apor- 
ten datos de carga hidráulica para lograr una 
mejor estimación en la zona de estudio. 

Es importante recalcar que la formulación 
del método se hizo para seleccionar en 
forma óptima tanto posiciones como tiempos 
de monitoreo, y la aplicación se hizo para 
la selección de posiciones de monitoreo 
únicamente. En un trabajo futuro se planea 
realizar una prueba similar para el diseño 
espacio-temporal de una red de monitoreo 
piezométrica. 

Finalmente es importante hacer notar que 
al aplicar este tipo de métodos a casos reales, 
además de utilizar criterios de optimización, es 
indispensable usar criterios geohidrológicos, 
que dependen del conocimiento del área de 
estudio. 
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Abstract 


BRISEÑO-RUIZ, J.V, HERRERA-ZAMARRÓN, G.S. € JÚNEZ-FERREIRA, H.E. 
Method for the optimal design ofnetworks to monitor groundwater levels. Water Technology 
and Sciences (in Spanish). Vol. 11, No. 4, October-December, 2011, pp. 77-96. 


It is essential to measure groundwater levels for long periods of time in order to estimate 
groundwater availability and its evolution. However, economic resources are a limitation 
to gathering enough data to estimate the evolution of these levels with certainty. Therefore, 
methods are required to design groundwater monitoring networks so as to effectively use the 
economic resources available for monitoring. The Herrera and Pinder method (Herrera, 1998) 
for the optimal design of groundwater quality sampling networks, applied by Herrera et al. 
(2001) and by Herrera and Pinder (2005), has several advantages: it allows for optimally 
selecting both monitoring positions and monitoring times; the information from the 
monitoring network enables obtaining estimates in space and time and; the method uses a flow 
and transport model that synthesizes knowledge about the characteristics of the aquifer, the 
groundwater flow dynamics and the groundwater quality problem, as well as knowledge about 
the laws of physics. This work tests, for the first time, a modification of the original method for 
the design of groundwater level monitoring networks. A synthetic case is presented in which 
the network is designed for a single time, assuming steady-state flow. The results indicate that 
the methodology proposed by Herrera and Pinder has good prospects for successful application 
to the optimal design of groundwater level monitoring networks. 


Keywords: groundwater, monitoring networks, optimal design. 
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